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Abstract 

We exactly calculate the fourth virial coefficient for hard spheres in even dimensions for D — 4, 6, 8, 10, 
and 12. 

Keywords: hard spheres, virial expansion. 

1 Introduction 

The virial series for the pressure 



P 



= P 



k B T 

a k=2 

of the system of hard spheres with diameter a in D dimensions specified by the two body pair potential 

has been studied for over 100 years. However, despite the long history of this problem, there are very few 
analytic results known. The second virial coefficient is 

B - aD * D/2 (3) 
B2 ~2r(l + D/2) [6) 

The third virial coefficient has been computed long ago for D = 2 £Q and D = 3 and in arbitrary 
dimension is compactly given as [3] 

4TY1 + R) t-x/3 
= ^l+D, i d<P (sin fl" (4) 







For low dimensions these are given in Tabled 

However, B4 has been calculated analytically in only two and three dimensions 

?± { 2 "^ + ^ for £> = 20||S1 (5) 

B * " I ^l-^ + ^arctan^ for D = 30E1I3 U 

and no further analytic results are known. An interesting discussion of the history of the calculation of B4 
in three dimensions is given in [Sj. 

We feel that for such a problem in pure geometry that tractable analytic results must exist. In this 
paper we make a modest step toward verifying this by extending the analytic evaluation of the fourth virial 
coefficient to even dimensions D < 12. 
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Table 1: The second and third virial coefficients as functions of dimension. 
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In the Mayer formalism 9, the fourth virial coefficient is 

b.— igi-|0-|n o» 

where each solid line represents the Mayer / function 

/(r)=cxp(-(7(r)/fc B T)-l (7) 
which for the hard sphere potential reduces to 

The second and third diagrams in this expansion have been evaluated in arbitrary dimension in [51 110|. In 
this paper we complete the computation in even dimensions for D < 12 by evaluating the first diagram in 
Section H 

For virial coefficients of order greater than four it is much more efficient to use the expansion of Ree and 
Hoover |11M12) where in addition to the / bonds we also have bonds /(r) = 1 + /(r) which are represented 
by dotted lines. In this notation every point is connected to every other point by either an / or an / bond 
and for the virial coefficient is given by 

where only the / bonds are shown. The integral in the first term is identical with the integral in the first 
term of Eq. and we compute the second term from the three Mayer diagrams of Eq. in Section [31 The 
results for B4 and for the two separate Ree-Hoover diagrams are given in Table |21 

2 Analytical Calculation of the Complete Star Diagram 

The complete star integral in the expansions of Eq. and Eq. [5] is by definition 



1 

V^oo V 



lim / d D rid D r2d D r3d D r4/(ri2)/(ri3)/(ri4)/(r23)/(r24)/(r34) 



d D r 1 d r2d r3/(r 1 2)/(r 1 3)/(r23)/(r 1 )/(r2)/(r 3 ) (10) 
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Table 2: Analytical results for the four point Ree- Hoover diagrams and B4 in even dimensions, with numerical 
values. 
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where ry = ri — r=. Specializing to the case of hard spheres we first note that /(r) = /(|r|), and then we 
treat the coordinates in Eq.^|as follows: ri is constrained to a unit ball centered on the origin due to / (fi), 
r 2 is integrated in the same ball with the additional condition that r%2 = |ria | = \P \ + r \ — 2r\r 2 cos 9 < 1. 
The r3 integral may be thought of as giving the overlapping hypervolume Vd of Z)-dimensional unit balls 
situated at the origin, rx, and r 2 , and this depends only on ri, r 2 , and the angle 9 between ri and r 2 . We 
define Vd as 

Vn(n,r 2 ,6) =-J -i D r 3 /(r 13 )/(r 23 )/(r 3 ) (11) 

where the negative sign ensures that Vd is positive. We may substitute this in to the expression for [x], and 
then write ri and r 2 in D-dimensional spherical polar coordinates 



)(2) 
L D-2 



d u v 1 d u v 2 f{r l )f{r 2 )f{r l2 )V D {r l ,r 2 ,9) 

= -J dr^- 1 J dn^J dr^- 1 J d9{sin 9)°- 2 J d£l\ 

xf(ri)f(r 2 )f(ri 2 )V D (ri,r 2 ,9) (12) 

where the angular integrals give Q.d-1 = J dSl D -i = 2tt d / 2 /T(D/2), which is valid for arbitrary D > 0, 
including non- integer D. The integrand in Eg. 1121 only depends on r\, r 2 , and 9, and we may integrate out 
the other angles to obtain 

^ = -n D _ifi D _ 2 J dnr?- 1 J dr 2 r®- x J d9{sin9) D - 2 

xf(ri)f(r 2 )f(r 12 )V D (ri,r 2 ,e) (13) 

We now change coordinates from (ri,r 2 ,9) to the coordinate system (R,ce,f3) which is illustrated in Fig. ^ 
The three points which were at the origin, ri , and r 2 are now circumscribed by a circle of radius R; 2a and 
2(3 are the angles subtended by r\ and r 2 from the center of the circle. 

r\ = 2R sin a 
r 2 = 2i?sin/3 

9 = n-a-P (14) 

Noting also that — f(ri)f(r 2 )f(ri 2 ) may only have the values of and 1, we may absorb the / functions in 
to the domain of integration and rewrite the triple integral of Eq. 1131 as 

da d/3 d J Ri? 2D ~ 1 [sinQsin/3sin(a + /3)] D_1 V D (J?,Q,/3) (15) 



3 



A 



3 2 



Figure 1: Change of variables. Figure 2: Domain of integration 

for £?4, where the shaded region 
satisfies the condition a > f3, 7. 



We are required to integrate over all R, a, and (3 such that the distance between any two points is less 
than one; denoting the sides opposite points A, B, and C on Fig. ^ as a, b, and c respectively we have 
a = 2i?sina, b = 2i?sin/3, and c = 2i?shi7. We choose to integrate over R before the angles a and /?, 
and we impose b, c < a < 1 which introduces an extra factor of 3. Therefore we must change the region of 
integration from 

0< b< a (16) 
0< c< a (17) 
0< a< 1 (18) 

to a domain for (R, a, (3). The boundary of this region is found by setting inequalities II fil and II 71 to equality, 
with solutions: 

a = (3 

a = 7r — [3 

a = 7 => /3 = 7r — 2a 

a = 7r-7^/3 = (19) 
This region is shown in Fig. [21 Finally, inequality 1181 restricts R to the domain 

< R < — — (20) 
2 sine* v ' 

The next step is to calculate the overlapping volume Vd of the three hyperspheres for even dimensions 
D = 2n, where n is an integer greater than or equal to 2. We do this in Appendix El and the result for 
R < 1 is denoted V£ and given in Eq. 0HI while the result for R > 1 is denoted Vp and given in Eq. I4*T1 
Throughout this paper we will adopt the convention that sums with index £ are always over {ot,/3, 7}. 

We will sketch the steps involved in the rest of the calculation, which were carried out using the computer 
algebra platform MAPLE version 8 for n — 2, 3, 4, 5, 6. The first step is to perform the integral over R, where 
we notice from inequality I2UI that we may only have R > 1 when a > Ignoring for the moment factors 
not involving R in Eq. we therefore need to calculate 



1/2 sin < ' 

6~ 



I 2 = I dR i? 4n -V< a < (21) 
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and 



rl/2s 



dR R 



All-\\r> 



V> a> 



5tt 



(22) 



h = J dR R in -W< + J 
One of the integrals we will need to do is 

= ^R 4n (arccos(i?sin£) - Rs\n£ yJ\-R? sin 2 ^ 

+ _L / dRR 4n sin£ Vl--R 2 sin 2 C 
2n 7 

and thus the only non- trivial parts of the integral are in the form J dR R 2l U m+ ^ with the shorthand notation 
U = 1 - R 2 sin 2 £. We substitute 

i? 2 = [1 — (1 — i? 2 sin 2 £)]/sin 2 £ everywhere, and this leaves us only to perform integrals such as 



(23) 



/ 



dR U j+ 2 



Vu 
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(24) 



We may replace 7 everywhere in this expression by (3, as we are integrating over a region that is symmetric 
in [3 and 7. After doing this and substituting in the limits of integration, we are left with terms involving 

[1 — sin 2 /3/(4sin 2 a)] 1 / 2 and arcsin[sin/3/(2 sin a)], along with functions of the form x~ m P(x) where m is 
an integer, P(x) is a polynomial, and x may be either sin a or sin/3. We now make the change of variables 
x = sin a, y = sin/3, and so the a and (3 integrals become 



da I dfi [sin a sin (3 sin(a + j3)] 



= I dx I ( 

J VT^J 
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The ambiguous sign is + when a < ^ and — when a > ^. Under the change of coordinates [1 — 
sin 2 /3/(4sin 2 a)] 1 / 2 becomes [4x 2 — y 2 ] 1 ^ 2 /(2x), and so naively it appears that we have to compute el- 
liptic integrals of the form J dz\Jt 2 — z 2 j\J\ — z 2 and then use identities for elliptic integrals to reduce the 
final result to a simple form. However, by splitting the final integrals over x and y for [X] in to the pieces 



dx x 



dyy 



2n-l T 2 



[x 2 (l-» a )] J [» a (l-x a )] 



2\1 n ~ 1— i 



(h or J 3 ) 



Jdxx 2 " 1 j dyy 2n 1 [y 2 (l - x 2 )] 3 [x 2 (l - y 2 )] 



2\1 " — 3 



(h or J 3 ) 



we avoid this. The first integral may be completed straightforwardly by first integrating over x, and then 
when one substitutes in the limits of integration this eliminates any elliptic integrals. The integral over y 
may now be performed without requiring any elliptic functions. Conversely, the second integral is completed 
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by first integrating over y and subsequently over x. Thus one can see that for any even dimension all integrals 
that must be performed are elementary. 

So we obtain the results of Table [5] by carrying out this procedure for D = 4,6,8,10, and 12 using 
MAPLE. 



3 Analytic derivation of the Ree-Hoover ring diagram 

The second and third Mayer diagrams have been found in terms of integrals over Bessel functions and 
hypergeometric functions by Luban and Baram The Mayer ring diagram is given by 

2 D + 4 r(l + D) [r(l + f)] 3 (l \-D_ 3 + D 3 + D \ 

U - r(i + ^)[r(§ + f)] 2 3 2 U' ' 2 ' 2 ' 2 'V 

while the second Mayer diagram is given in as 

= -2 D+1 D 3 [T(D/2)} 2 f dy y [g D/2 (y)} 2 

Jo 

where 

g v (y)= / dx x [J v (x)\ J v -i(xy) 



Other expressions for these diagrams were given in but MAPLE was straightforwardly able to evaluate 

these expressions for dimensions one through twelve, and these are listed Table |3J The second diagram of 
the Ree-Hoover expansion was then obtained from the equation 

! = ^ + 20 + n (25) 



Table 3: and | in dimensions up to twelve. 
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A Calculation of Vd in even dimensions 



For D — 2 the calculation of Vd is straightforward as the r3 integral is confined to the same plane as the 
three circle centers and trivially gives the area of intersection of the three circles of radius r = 1 . Although 
we will neglect the details for the moment, the area depends on R, a, and (3 and we will denote this as 
A(r = l,R,a,0). 

As shown in Fig. [3 for D = 3, we define the perpendicular distance from the plane of the circle centers 
to be h; h = is the original plane in which we see three intersecting circles of radius one. As we increase 
h we see overlapping circles with the same center but decreasing radius, with the radii of these circles given 
by r = \/\ — h 2 (see Figs. 0]and|3J). The total volume of intersection may be obtained by integrating the 
overlapping area of the three circles with respect to h, from h = to the value h max where the intersection 
of the three circles is reduced to a single point. There is an additional factor of two because we need to 
integrate both above and below the plane. 



V 3 (R,a,P) = 2 / dhA(r,R,a,/3)\ r = VI= jp 
Jo 



= 2 



dr 



A(r,R,a,f3) 



(26) 




Figure 3: Intersecting three dimensional spheres. 



For the general case of D > 3, we choose a spherical polar coordinate system for the D — 2 dimensional 
subspace with radial coordinate h; there is an extra factor of h D ~ 3 compared to three dimensional case, and 
we now need to perform an additional trivial angular integration, ie 



V D (R,a,/3) = 



dhh 



D-3 



dn D _ 3 A(r, R, a, 0) \ r= ^W 



9. 



D-3, 



dhh u - 6 A{r,R,a,P)\ r=VI=w 
drr{l~r 2 ) {D - i)/2 A{r, R,a,/3) 



(27) 



Although Eq. |27| breaks down for D < 2, note that it is correct for D ~ 3 as fio = 2. 

As the area of intersection may be bounded by arcs from either two or three circles, we will need to deal 
with two different expressions for the area, and three different cases for Vd- In Fig. [3] we have R < 1 and 
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Figure 4: Circles produced by 
cutting intersecting spheres at 
h = 0, h = 1/V5, and h = 
V 1 — R 2 = v3/2, where circle 
centers form an acute angled tri- 
angle. 



Figure 5: Circles produced by 
cutting the intersecting spheres 
at h = 0, h = 1/V%, and h = 
y/1 - R 2 = V3/2. In this case 
the circle centers form a trian- 
gle with an obtuse angle and the 
area of intersection will be non- 
zero even for h > \/l — R 2 . 



a < \ in which case the area is bounded by three arcs over the full range of r, while in Fig. [3] we have R < 1 
and a > J f° r which the area is bounded by three arcs for R <r < 1, but by two arcs for i?sina < r < R. 
When R > 1 the area can only have two arcs for its boundary and this is the third case. 

Provided the overlapping region is bounded by arcs from three circles, the area of overlap may be 
calculated by an inclusion-exclusion method following the method used for the overlapping volume of 3 
spheres in |13| . The area of intersection of three circles of radius r is given by 

A^(r,R 7 a,(3) 

= (Area of A ABC) — (Area of sector within circle Aj) 

Ai = {A,B,C} 

+ - (Area of intersection of circles A; and Aj) 

{Ai,Aj} 

= ii? 2 (sin 2a + sin 2/3 + sin 27) - ^r 2 (a + /? + 7) 



i?sin£ 



r 

1 



^ ( r 2 arccos ( ) — R sin £ \j r 2 — R 2 sin 2 £ 



= i? 2 ^sin£cos£ — -irr 2 

+ [r 2 arccos ) - i?sin£ yj ' r 2 - R 2 sin 2 £j (28) 

The different contributions to Eq. EHlare shown explicitly in Fig. El and one should also note that although 
expressions are written in terms of a, (3, and 7 for the sake of simplicity, these variables are not independent 
and a + /3 + 7 = tt. The overlap of two circles of radius r separated by a distance i?sina is given by 



A^> (r, R, a) = 2 ^r 2 arccos ^ — S * n - ^ - R sin a V V 2 - R 2 sin 2 a ^ 



(29) 



We may obtain this expression from Eq. |2HI by taking the limit r — * R + when a > ^ . Thus the three 
expressions for Vd are 

r " / ' 1 drr(l-r 2 )( D - 4 )/ 2 ^ 3 ) «<J (30) 



" drr (l_ r 2)(0-4)/ 2j4 (3) 



^D-3 Jr 



8 



Figure 6: The different kinds of contributions to the inclusion-exclusion formula for the area of intersection 
of three circles: (a) Area of AABC (b) Area of sector A (c) Half of the area of intersection of circles A and 
B. 



V D 



D-3 



+ T drr(l^ 2 )< D - 4 )/M< 2 ) a>^R<l 

JRsina ^ 

f drr{l-r 2 y D -^' 2 A^ a > J R > 1 



(31) 
(32) 



We will now specialize to the case where D — 2n, where n is an integer greater than or equal to 2. By 
defining 



h(r,0 



f 

Jo 



dr r 3 (l - r 2 ) r 



^fisin^ i?sin£ ^ ^i?sin^ 



and integrating out the remaining pieces, we have for a < n/2, R < 1: 



V D 



Or 



dr r(l — r 2 ) n 2 



R 2 sin g cos g — i-; 



(l-R 2 



2(n- 1) 
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2n 
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+ 2(M 1 >0--Ti(fl.O) 



while for a > 7r/2, i? < 1: 

f2_D_3 



2(n-l) 





-i? 2 sing cos £ + 



(l-fi 2 )" _ (l-J? 2 )"- 1 ^ 7T 

2n 2(n - 1) J 2 



+ X]( 7l ( 1 '^ - h ( R ^0) + 2(/i(ii,a) - 7i(iisina,a)) 



(1 - 7? 2 ) r 



(l-i? 2 ) n (l-i? 2 



2 { n-— R2 ^^ COSi+ \ 2n 2(n - 1) ^2 

-^7i(l,0+/i(-R,a) -h(R,(3) -h(R,j) - 2h(Rsma,a)) 



and finally for a > ^ i? > 1: 



Vb 



n 



D-3 



2(7i(l,a) - h(Rsma,a)) 
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Integrating I\ once by parts 



h(r,0 = 



(1-r 2 



+ 



1 



2n 2(n - 1) 2n(n - 1) 

' R sin £ \ i? sin £ 



/i?sin£\ 



Y^tlL( n - 2 \ / drr 2 '' +1 _i_£ V^ 2 - i? 2 sin 2 \ 
, =0 J + 2 V J / -A. 



(l-r 2 r 
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1 
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' R sin £ \ i? sin £ 
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E 
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{-Vf_(n-2 
7+2" 



/■r 

/ dr r((r 2 - R 2 sin 2 £) + i? 2 sin 2 £) j .Rsin£ y/r 2 - R 2 sin 2 £ 
(1-r 2 )"- 2 (1-r 2 )™- 1 



+ 
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' _R sin £ \ i? sin £ 



^ fisin^ _ Rsm£ ^ .Bsin^ y 



J +2 



(7?sin£) 



2j-2fe+l 



x [ drr(r 2 ~R 2 sin 2 £) fe+ * 
Jo 

(1-r 2 )™- 2 (1-r 2 )™" 1 



+ 



2n 2(n - 1) 2n(n - 1) 

' i? sin £ \ i? sin £ 



/flsmf \ 

■(— )- 



n-2 j 



J + 2 v J 

(r 2 -i? sin 2 C) fc+ * 



j=0 k=0 
- 2 Tf2-2^k+ 



2k + 3 



where we have used 



/' 

Jo 



drr 6 (l-r 2 ) 



2\n-2 



drr(l-(l-r 2 ))(l-r 2 ) 

(1-r 2 )"' 1 (1-r 2 )" 
_____ + _ . |_ 



2\n-2 



2n 



2n(n- 1) 



n-2 

IE 



(-lp /n-2 



.2.7+4 



(33) 



(34) 
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and 



d 
dr 



arccos 



i?sin^\ i?sin£ 



i?sin£ 



2i?sin£ 



'1- 



i?sin£ 



Substituting in the limits of integration, r = 1 



|arccos(i?sin£) — ijsinf \/l — R 2 sin 2 £ 

2n(n - 1) L 



n-2 3 
3=0 fc=0 



2 7 - 2fc +l (l-^ 2 sin 2 g)*+* 
2fc + 3 



Substituting r = R 



(1 - fl 2 )* 



(1 - K 2 ) r ' 



■ + ■ 



2n 2(n - 1) 2n(n - 1) 

j^arccos(sing) — sin£-\/ cos 2 £ 



ra — 2 3 



3=0 fc=0 
sgn 



j + 2 v j ' W 



(l-i? 2 )"- 2 (l-R 2 )™- 1 



2/c + 3 



(H 



2n 



2(n - 1) 2n(n - 1) 



X |^ - £ - sin £ cos £ 

n-2 j 



-E£^(V)0^ +4(sin ^~ 2fe+1 

3=0 fc = 

= sgn(|-?) /<(/?,£) 

which defines I K (R, £) . Substituting r = R sin £ 

h(R S m£,£) = 



(cos £) 



2fe+3 



2fc + 3 



(35) 



(36) 



(37) 



(38) 



Substituting Ii(R,^) — sgn(^ — in to Vd, we notice that the two cases of a < n/2 and 

a > tt/2 with R < 1 reduce to the same form, which we will denote as Vp . Similarly denotes the volume 
of overlap when R > 1. 



Or 



+5^(/i(i,o-/i < (fl,o) «<i 



(i - i? 2 )" (l-i? 2 )" 



2n 



2(n- 1) 



*2 



27i(l,a) R>1 

Finally, by substituting the expressions for I\ and using a + /3 + 7 = 7rwe obtain for R < 1 : 

. ^— — -y ^arccos(i?sin£) — i?sin£ \/ 1 — R 2 sin 2 ^ 



Or 



(39) 



Hi 



11 



and for R > 1: 



n-2 j 

EEE ~J+2 



(-iy fn-2\ /.A /D „.„ rt 2,- a+ i(l-fl 2 sin 2 C) t+ ' 



5 j=0 fc=0 

n-2 j 

+EEE 



(i?sin£) 



2fc + 3 



(-1)3 / n - 2^ ^ D 2,+4,„ : „ c ,2,-2 fc + l(cOs0 2fc+ ; 

j + 2 I i 



5 j=0 fc=0 

(1- i? 2 ) n 1 

2n(n - 1) + 2n(n - 1) 



r,27+4/ • ^\2j-2fc + l 1 

fc y* (sm?) -2* + 3- 



- — — ■ ^arccos(_Rsina) — i? sin a \/l — R 2 sin 2 



n(n — 1) 

n-2 j 



; EE 



J + 2 



j=0 fc=0 

(1 - fl 2 sin 2 a) fc +i 
2fc + 3 



(i?sina 



2j-2fe+l 



(40) 



(41) 
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